from IPython.display import display
from IPython.display import HTML
import IPython.core.display as di # Example: di.display_html('<h3>%s:</h3>' % str, raw=True)
# This line will hide code by default when the notebook is exported as HTML
di.display_html('<script>jQuery(function() {if (jQuery("body.notebook_app").length == 0) { jQuery(".input_area").toggle(); jQuery(".prompt").toggle();}});</script>', raw=True)
# This line will add a button to toggle visibility of code blocks, for use with the HTML export version
di.display_html('''<button onclick="jQuery('.input_area').toggle(); jQuery('.prompt').toggle();">Toggle code</button>''', raw=True)
%reset
#pylab qt -> open windows, needs display
%pylab qt
%load_ext autoreload
%autoreload 2
#from pySurf.fit_cylinder import *
#from pySurf.points import *
from pySurf.points import plot_points,matrix_to_points2
from pySurf.points import subtract_points,level_points
from pySurf.points import points_find_hull,crop_points
from pySurf.psd2d import calculatePSD
from pySurf.points import points_autoresample
#from calibrate_align import *
from pyProfile.profile import polyfit_profile
from plotting.multiplots import compare_images
from dataIO.fn_add_subfix import fn_add_subfix
from dataIO.span import span
#from PSDanalysis import *
from pySurf.instrumentReader import matrix4D_reader
import os
from pySurf.instrumentReader import bin_reader
from plotting.add_clickable_markers import add_clickable_markers
from plotting.backends import maximize
from pySurf.instrumentReader import points_reader
from pySurf.affine2D import plot_transform,find_rototrans
from pySurf.distanceTable import distanceTable
from pySurf.data2D_class import Data2D
from pySurf.data2D import data_from_txt
from dataIO.outliers import remove_outliers
%qtconsole
pwd
from pySurf.psd2d import calculatePSD
from utilities.imaging import fitting as fit
correct_misalignment=fit.fitConeMisalign
#2019/03/19 C1S04_PZT_stress_fit from C1S05_stress_compensated_WFS_fit
# fit PZT stress on first measuremnt of sample received from PSU,
# gold coating applied at PSU.
outfolder=r'results\C1S04_PZT_stress_fit' #set output
WFS data are utilized because the deformation due to PZT stress is likely to be out of the 4D interferomenter dynamic range. Data are presented below:
wf=r'G:\My Drive\POOL\PSD\analysis_PSD\results\C1S04_cut\180206_C1S04_GentexCut_Meas3_matrix.dat'
p1=Data2D(file=wf,reader=data_from_txt,delimiter=' ',units=['mm','mm','$\mu$m'],
name='C1S04 substrate WFS data',addaxis=1,center=(0,0))
plt.figure(1)
plt.clf()
p1.level().plot()
display(plt.gcf())
wf=r'G:\My Drive\POOL\PSD\analysis_PSD\results\C1S04_PZT\190314_06_C1S04.dat'
p2=Data2D(file=wf,reader=data_from_txt,delimiter=' ',units=['mm','mm','$\mu$m'],
name='C1S04 PZT WFS',addaxis=True,center=(0,0))
p2=p2.remove_nan_frame()
plt.figure(1)
plt.clf()
p2.plot()
display(plt.gcf())
Data are aligned by interactively selecting recognizeable features and dimples. An affine transformation is employed to correct potential difference of horizontal scale.
Picture below shows data before and after PZT with first 2 2D Legendre terms removed.
plt.close('all')
plt.figure()
maximize()
ax1=plt.subplot(121)
p1.level(2).plot()
plt.clim([-1,1])
ax2=plt.subplot(122,sharex=ax1,sharey=ax1)
p2.level(2).plot()
plt.clim([-1,1])
plt.suptitle('2nd order 2D legendre removed')
display(plt.gcf())
plt.savefig(os.path.join(outfolder,'data_comparison_2leg.png'))
Here same data with sag removed by line, as used for feature identification and alignment:
plt.figure()
maximize()
ax1=plt.subplot(121)
p1.level(2,byline=True).plot()
plt.clim(*remove_outliers(p1.level(2,byline=True).data,span=True,
nsigma=2,itmax=3))
ax2=plt.subplot(122,sharex=ax1,sharey=ax1)
p2.level(2,byline=True).plot()
plt.clim(*remove_outliers(p2.level(2,byline=True).data,span=True,
nsigma=2,itmax=3))
plt.suptitle('Y sag removed by line - natural color scale')
display(plt.gcf())
plt.savefig(os.path.join(outfolder,'data_comparison_sag.png'))
from plotting.add_clickable_markers import add_clickable_markers2
add_clickable_markers2(ax=ax1)
add_clickable_markers2(ax=ax2)
These values are the markers positions:
print(ax1.markers)
print(ax2.markers)
m1=np.array([[24.58277715001921, 21.497178162325042], [19.29575731191258, 28.96120616906383], [-3.251827291777502, 28.805705585590104], [18.051752644122786, -11.624446117578294], [10.121222886962826, -17.377967706106105], [-46.32548891399922, 32.84872075590694], [-38.39495915683926, 41.401252846961796], [-43.059976661051, -26.086000380634687]])
m2=np.array([[-20.15881919858859, -20.478152809990057], [-14.408057148043014, -27.9386008755627], [7.507009044576591, -27.62774887283051], [-15.340613156239613, 12.627585480988529], [-6.947609082470393, 18.689199534266308], [45.43095337790416, -30.73626890015244], [39.05848732189423, -40.06182898211824], [41.70072934511788, 26.305073601205038]])
To verify the scale of the images, the three tables below report the distance of each marker from all the others on same surface data.
These are reported respectively for each of the data. The third table is the difference in corresponding distances.
Picture below shows the position of the fiducials used for alignment.
from pySurf.distanceTable import distanceTable
#import print_options
#print_options.set_float_precision(2)
print(distanceTable(m1))
print('---------------------------')
print(distanceTable(m2))
print('---------------------------')
print(distanceTable(m1)-distanceTable(m2))
plt.close('all')
plt.figure()
maximize()
plt.suptitle('2 legendre removed by line for alignment')
ax1=plt.subplot(121)
p1.level(2,byline=True).plot()
plt.clim([-1,1])
ax2=plt.subplot(122,sharex=ax1,sharey=ax1)
p2.level(2,byline=True).plot()
plt.clim([-1,1])
for a,m in zip([ax1,ax2],[m1,m2]):
plt.sca(a)
for mm in m:
plot(mm[0],mm[1],'rx')
display(plt.gcf())
plt.savefig(os.path.join(outfolder,'data_markers.png'))
import logging
try:
logger
except NameError:
# create logger
logger = logging.getLogger()
logger.setLevel(logging.DEBUG)
# create console handler and set level to debug
ch = logging.StreamHandler()
ch.setLevel(logging.INFO)
# create formatter
formatter = logging.Formatter('%(asctime)s - %(name)s - %(levelname)s - %(message)s')
# add formatter to ch
ch.setFormatter(formatter)
# add ch to logger
logger.addHandler(ch)
Two transformations are chosen: rototranslation and affine transformation, that transform second set of markers (on sample after PZT) to the other.
In an ideal case, transforming markers should bring to exact same positions on first untransformed data. Errors on markers position, respectively with rototranslation and affine tranformation, are reported below. Errors include optical (metrology) effects and the inaccuracy in the selection of exact marker location.
Affine transformation is clearly more accurate and it will be used to align the two data (the second one is transformed).
from pySurf.affine2D import find_affine
trans=find_rototrans(m2,m1)
transa=find_affine(m2,m1)
p2t=p2.apply_transform(transa)
mm2=transa(m2)
Table of distances with transformed data show that the transformation effectively compensated the scale errors:
from pySurf.distanceTable import distanceTable
#import print_options
#print_options.set_float_precision(2)
print(distanceTable(np.array(m1)))
print('---------------------------')
print(distanceTable(np.array(mm2)))
print('---------------------------')
print(distanceTable(np.array(m1))-distanceTable(np.array(mm2)))
The resulting net PZT effect (PZT shape - substrate shape) is shown below as surface and central axial profile.
def diffplot(p1t,p4):
"""plots and return difference"""
plt.figure()
maximize()
ax1=plt.subplot(131)
p1t.level((1,1)).plot()
plt.title('PZT')
plt.clim([-3,3])
ax2=plt.subplot(132,sharex=ax1,sharey=ax1)
p4.level((1,1)).plot()
plt.title('substrate')
plt.clim([-3,3])
ax3=plt.subplot(133,sharex=ax1,sharey=ax1)
diff=(p1t-p4).level((1,1))
diff.name='PZT net effect (PZT-substrate)'
diff.plot()
plt.clim([-1,1])
display(plt.gcf())
return diff
diffPZT=diffplot(p2t,p1)
p1.save(os.path.join(outfolder,'wfs_substrate_aligned.dat'))
p2t.save(os.path.join(outfolder,'wfs_PZT_aligned.dat'))
diffPZT.save(os.path.join(outfolder,'figure_change_aligned.dat'))
plt.savefig(os.path.join(outfolder,'figure_change_aligned.png'))
from pyProfile.profile import level_profile, subtract_profiles
from pySurf.points import extract_profile
## Extract axial profile
x1,y1 = level_profile(*(extract_profile(p1.topoints(),(0,-45),(0,45))))
x2,y2 = level_profile(*(extract_profile(p2t.topoints(),(0,-45),(0,45))))
xd,yd = level_profile(*subtract_profiles(x2,y2,x1,y1))
plt.figure()
plt.plot(x1,y1,label='Substrate')
plt.plot(x2,y2,label='PZT')
plt.plot(xd,yd,label='Profile change')
plt.xlabel(p1.units[0])
plt.ylabel(p1.units[2])
plt.grid(1)
plt.ylim([-0.8,0.8])
plt.legend(loc=0)
plt.title('central axial profiles')
display(plt.gcf())
plt.savefig(os.path.join(outfolder,'data_central_profile.png'))
Fit of stress is calculated starting fomr Vladimir Kradinov's simulation 2018/06/19 corresponding to a coating with uniform integrated stress 300 MPa um.
This is shown below alone and in comparison with data.
from pySurf.instrumentReader import points_reader, fitsWFS_reader
from pySurf.affine2D import plot_transform,find_rototrans
from pySurf.distanceTable import distanceTable
from pySurf.data2D_class import Data2D
from pySurf.data2D import resample_data
FEAf=r'G:\My Drive\Shared by Vincenzo\Metrology logs and data\Simulations\Coating_stress\VK20180619\coat_300_MPa_um_ur.fits' #300 MPa*um load
plt.figure(2)
plt.clf()
#manually adjust x and y scales to 4 in, adjsut units and signs.
fdata,xf,yf=fitsWFS_reader(FEAf) #can accept ypix, ytox,..
#xf,yf,fdata=yf,xf,fdata.T
fdata=fdata-correct_misalignment(fdata)[0]
fdata=-fdata*1000000. #m to um
#fitsWFS_reader makes a mysterious change of sign, that is good with WK's
#data that are bump negative conversion.
#scale and center axes (100 mm square image size)
nx=fdata.shape[0]
scale=101.6/nx
xf=(xf-span(xf).sum()/2)*scale
yf=(yf-span(yf).sum()/2)*scale
df=Data2D(fdata,xf,yf,units=['mm','mm','$\mu$m'],name='VK Simulation (cone)')
#rule is: corners have same sign as delta Radius
plt.clf()
df.plot()
display(plt.gcf())
df.save(fn_add_subfix(FEAf,'','.dat',strip=True,pre=outfolder+os.path.sep))
plt.savefig(fn_add_subfix(FEAf,'','.png',strip=True,pre=outfolder+os.path.sep))
fn_add_subfix(FEAf,'','.dat',strip=True,pre=outfolder+os.path.sep)
p1=diffPZT
diff=p1-df
plt.figure()
maximize()
ax1=plt.subplot(131)
p1.plot()
ax2=plt.subplot(132,sharex=ax1,sharey=ax1)
df.plot()
ax3=plt.subplot(133,sharex=ax1,sharey=ax1)
diff.plot()
display(plt.gcf())
plt.savefig(os.path.join(outfolder,'starting_data.png'))
outfolder
from stress_compensation import comp2, scale_fit2, trioptimize, fom, fom_sl
import pandas as pd
The best fit for stress is determined by raster scan of a scaling value for the FEA data and chosing the value that minimizes the rms of surface difference with PZT net effect. The raster scan is executed twice, coarse (no output) and finer (shown below), to determine the scale factor with two decimal digits.
The fit is repeated for different azimuthal apertures. For each of them it is shown:
Each plot reports values for the ROI (small box) and for the entire surface (larger box) (the two are identical for the crop image at point 2)
Used functions (for who is interested in details).
print(help(trioptimize))
print('---------------------------')
print(help(scale_fit2))
print('---------------------------')
print(help(fom))
print('---------------------------')
print(help(fom_sl))
allres=pd.DataFrame()
os.makedirs(os.path.join(outfolder,'fits'),exist_ok=True)
#ROUGH FIT - set dis to True in call to show results
scales=trioptimize(p1,df,fom=fom,testvals=[np.arange(0.1,0.7,0.1),
np.arange(0.1,0.7,0.1),
np.arange(0.4,1.0,0.1),
np.arange(0.4,1.0,0.1),
np.arange(0.6,1.1,0.1)],
rois=[[[-10,10],[-45,45]],[[-20,20],[-45,45]],[[-30,30],[-45,45]],
[[-37,37],[-45,45]],None],dis =False,
outfile=os.path.join(outfolder,'fits','wfs_test_roi'))
#allres = allres.append(scales)
scales
#FINE TUNING
allres=trioptimize(p1,df,fom=fom,testvals=[np.arange(0.45,0.55,0.01),
np.arange(0.35,0.45,0.01),
np.arange(0.7,0.8,0.01),
np.arange(0.8,0.9,0.01),
np.arange(0.95,1.05,0.01)],
rois=[[[-10,10],[-45,45]],[[-20,20],[-45,45]],[[-30,30],[-45,45]],
[[-37,37],[-45,45]],None],
outfile=os.path.join(outfolder,'fits','wfs_test_roi'))
allres[['tbest','rms_roi','roi']]
allres[['tbest','rms_roi','roi']].to_csv(os.path.join(outfolder,'all_stats.txt'))
plt.close('all')
comp2(p1,df*1.02)
plt.suptitle('best factor = 1.02')
display(plt.gcf())
plt.savefig(os.path.join(outfolder,'best_fit_VK.png'))
The exact same procedure is repeated on a different set of data by Vanessa Marquez, representing a PZT film on back of substrate, having an integrated stress of 150 MPa um. It is then expected to find values that are about double of previous ones.
from pySurf.instrumentReader import FEAreader
FEAf2=r'G:\My Drive\Shared by Vincenzo\Metrology logs and data\Simulations\CoatingStress_Analysis\220mm_CoatingStress\Primary\220mmROC_Case1_MirrorSurfaceOnly.csv'
plt.figure(2)
plt.clf()
#manually adjust x and y scales to 4 in, adjsut units and signs.
fdata,xf,yf=points_autoresample(FEAreader(FEAf2)) #can accept ypix, ytox,..
#xf,yf,fdata=yf,xf,fdata.T
fdata=fdata-correct_misalignment(fdata)[0]
fdata=fdata*1000. #m to um
#fitsWFS_reader makes a mysterious change of sign, that is good with WK's
#data that are bump negative conversion.
#scale and center axes (100 mm square image size)
nx=fdata.shape[0]
scale=1.016
xf=(xf-span(xf).sum()/2)*scale
yf=(yf-span(yf).sum()/2)*scale
df=Data2D(fdata,xf,yf,units=['mm','mm','$\mu$m'],name='VM simulation 220mmROC_Case1_MirrorSurfaceOnly')
#rule is: corners have same sign as delta Radius
plt.clf()
df.plot()
display(plt.gcf())
df.save(fn_add_subfix(FEAf2,'','.dat',strip=True,pre=outfolder+os.path.sep))
plt.savefig(fn_add_subfix(FEAf2,'','.png',strip=True,pre=outfolder+os.path.sep))
os.makedirs(os.path.join(outfolder,'fits2'),exist_ok=True)
#ROUGH FIT
scales=trioptimize(p1,df,fom=fom,testvals=[np.arange(0.8,1.3,0.1),
np.arange(0.5,1.1,0.1),
np.arange(0.8,1.4,0.1),
np.arange(1.3,1.9,0.1),
np.arange(1.8,2.4,0.1)],
rois=[[[-10,10],[-45,45]],[[-20,20],[-45,45]],[[-30,30],[-45,45]],
[[-37,37],[-45,45]],None],dis=False,
outfile=os.path.join(outfolder,'fits2','wfs_test_roi'))
#allres = allres.append(scales)
scales
#FINE TUNING
scales=trioptimize(p1,df,fom=fom,testvals=[np.arange(0.88,0.92,0.01),
np.arange(0.64,0.71,0.01),
np.arange(1.09,1.14,0.01),
np.arange(1.38,1.42,0.01),
np.arange(1.92,2.0,0.01)],
rois=[[[-10,10],[-45,45]],[[-20,20],[-45,45]],[[-30,30],[-45,45]],
[[-37,37],[-45,45]],None],
outfile=os.path.join(outfolder,'fits2','wfs_test_roi'))
allres=allres.append(scales)
This is a summary of results (first half from VK's simulation, 300 MPa um, second half for VM's, 150 MPa um)
allres[['tbest','rms_roi','roi']]
allres[['tbest','rms_roi','roi']].to_csv(os.path.join(outfolder,'all_stats_VKVM.txt'))
bf=1.94
plt.close('all')
comp2(p1,df*bf)
plt.suptitle('best factor = %F.2'%bf)
display(plt.gcf())
plt.savefig(os.path.join(outfolder,'best_fit_VM.png'))
from pySurf.instrumentReader import points_reader, fitsWFS_reader, FEAreader
from pySurf.affine2D import plot_transform,find_rototrans
from pySurf.distanceTable import distanceTable
from pySurf.data2D_class import Data2D
from pySurf.data2D import resample_data
FEAf=r'G:\My Drive\Shared by Vincenzo\Metrology logs and data\Simulations\Coating_stress\VK20180619\coat_300_MPa_um_ur.fits' #300 MPa*um load
plt.figure(2)
plt.clf()
#manually adjust x and y scales to 4 in, adjsut units and signs.
fdata,xf,yf=fitsWFS_reader(FEAf) #can accept ypix, ytox,..
#xf,yf,fdata=yf,xf,fdata.T
fdata=fdata-correct_misalignment(fdata)[0]
fdata=-fdata*1000000. #m to um
#fitsWFS_reader makes a mysterious change of sign, that is good with WK's
#data that are bump negative conversion.
#scale and center axes (100 mm square image size)
nx=fdata.shape[0]
scale=101.6/nx
xf=(xf-span(xf).sum()/2)*scale
yf=(yf-span(yf).sum()/2)*scale
df=Data2D(fdata,xf,yf,units=['mm','mm','$\mu$m'],name='VK Simulation (cone)')
#rule is: corners have same sign as delta Radius
FEAf2=r'G:\My Drive\Shared by Vincenzo\Metrology logs and data\Simulations\CoatingStress_Analysis\220mm_CoatingStress\Primary\220mmROC_Case1_MirrorSurfaceOnly.csv'
plt.figure(2)
plt.clf()
#manually adjust x and y scales to 4 in, adjsut units and signs.
fdata,xf,yf=points_autoresample(FEAreader(FEAf2)) #can accept ypix, ytox,..
#xf,yf,fdata=yf,xf,fdata.T
fdata=fdata-correct_misalignment(fdata)[0]
fdata=fdata*1000. #m to um
#fitsWFS_reader makes a mysterious change of sign, that is good with WK's
#data that are bump negative conversion.
#scale and center axes (100 mm square image size)
nx=fdata.shape[0]
scale=1.016
xf=(xf-span(xf).sum()/2)*scale
yf=(yf-span(yf).sum()/2)*scale
df2=Data2D(fdata,xf,yf,units=['mm','mm','$\mu$m'],name='VM simulation 220mmROC_Case1_MirrorSurfaceOnly')
#rule is: corners have same sign as delta Radius
plt.close('all')
df=df.crop([-45,45],[-45,45]).level()
df.plot()
display(plt.gcf())
plt.figure()
df2=df2.crop([-45,45],[-45,45]).level()
df2.plot()
display(plt.gcf())
plt.figure()
maximize()
ax1=plt.subplot(131)
df.plot(title='VK')
plt.subplot(132,sharex=ax1,sharey=ax1)
(df2).plot(title='VM')
plt.subplot(133,sharex=ax1,sharey=ax1)
(df-df2*2).plot(title='VK-VM*2')
display(plt.gcf())
plt.savefig(os.path.join(outfolder,'FEAdiff_crop.png'))
plt.close('all')